## ANALYSIS REPLICATION FILE
# Paper: Coethnicity and Corruption: Field Experimental Evidence from Public Officials in Malawi
# Authors: Brigitte Seim and Amanda Lea Robinson
# Journal: Journal of Experimental Political Science

## FILE INFO
# This R code performs analyses reported in the appendix.
# Place this file in the same folder as the ESCOMdata_SR.csv data file and set that file as your working directory.

###########################################################

## Data Input and Preparation ##

rm(list = ls())
options(stringsAsFactors = FALSE)
data2 <- read.csv("ESCOMdata_SR.csv")
library(rstan)
library(brms)
library(coda)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
library(brms)
data2$outcome2 <- "Normal"
data2$outcome2[data2$outcome==0] <- "Bribery"
data2$outcome2[data2$outcome==2] <- "Expedited"
data2$outcome2 <- as.factor(data2$outcome2)
data2$outcome2 <- relevel(data2$outcome2, "Normal")
data2$SES <- data2$ses
data2$Power <- data2$power
data2$Coethnic <- data2$coreg
data2$NumberOfficials <- data2$numofficials
data2$OtherCustomers <- data2$other_customers

## Main Model ##

bayesescom3<-brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 10)), seed = 614)
summary(bayesescom3)
plot(bayesescom3,pars=parnames(bayesescom3)[3:5])
plot(bayesescom3,pars=parnames(bayesescom3)[8:10])
x3<-as.mcmc(bayesescom3)
HPDinterval(x3,0.9)
xall3<-do.call(rbind,x3)
HPDinterval(mcmc(xall3),0.9)

## Prior Sensitivity Analysis ##

# Models #
bayesescom4 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 9.8)), seed = 614)
bayesescom5 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 9.6)), seed = 614)
bayesescom6 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 9.4)), seed = 614)
bayesescom7 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 9.2)), seed = 614)
bayesescom8 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 9)), seed = 614)
bayesescom9 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 8.8)), seed = 614)
bayesescom10 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 8.6)), seed = 614)
bayesescom11 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 8.40000000000001)), seed = 614)
bayesescom12 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 8.20000000000001)), seed = 614)
bayesescom13 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 8.00000000000001)), seed = 614)
bayesescom14 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 7.80000000000001)), seed = 614)
bayesescom15 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 7.60000000000001)), seed = 614)
bayesescom16 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 7.40000000000001)), seed = 614)
bayesescom17 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 7.20000000000001)), seed = 614)
bayesescom18 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 7.00000000000001)), seed = 614)
bayesescom19 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 6.80000000000001)), seed = 614)
bayesescom20 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 6.60000000000001)), seed = 614)
bayesescom21 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 6.40000000000001)), seed = 614)
bayesescom22 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 6.20000000000001)), seed = 614)
bayesescom23 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 6.00000000000001)), seed = 614)
bayesescom24 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 5.80000000000001)), seed = 614)
bayesescom25 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 5.60000000000002)), seed = 614)
bayesescom26 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 5.40000000000002)), seed = 614)
bayesescom27 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 5.20000000000002)), seed = 614)
bayesescom28 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 5.00000000000002)), seed = 614)
bayesescom29 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 4.80000000000002)), seed = 614)
bayesescom30 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 4.60000000000002)), seed = 614)
bayesescom31 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 4.40000000000002)), seed = 614)
bayesescom32 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 4.20000000000002)), seed = 614)
bayesescom33 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 4.00000000000002)), seed = 614)
bayesescom34 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 3.80000000000002)), seed = 614)
bayesescom35 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 3.60000000000002)), seed = 614)
bayesescom36 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 3.40000000000002)), seed = 614)
bayesescom37 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 3.20000000000002)), seed = 614)
bayesescom38 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 3.00000000000002)), seed = 614)
bayesescom39 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 2.80000000000003)), seed = 614)
bayesescom40 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 2.60000000000003)), seed = 614)
bayesescom41 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 2.40000000000003)), seed = 614)
bayesescom42 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 2.20000000000003)), seed = 614)
bayesescom43 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 2.00000000000003)), seed = 614)
bayesescom44 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 1.80000000000003)), seed = 614)
bayesescom45 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 1.60000000000003)), seed = 614)
bayesescom46 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 1.40000000000003)), seed = 614)
bayesescom47 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 1.20000000000003)), seed = 614)
bayesescom48 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 1.00000000000003)), seed = 614)
bayesescom49 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 0.800000000000029)), seed = 614)
bayesescom50 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 0.60000000000003)), seed = 614)
bayesescom51 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 0.400000000000031)), seed = 614)
bayesescom52 <- brm(outcome2 ~ Coethnic + SES + Power + NumberOfficials + OtherCustomers, data = data2, family = categorical(link = "logit"),prior =  prior(normal(0, 0.200000000000029)), seed = 614)

# 90% Intervals #

x3<-as.mcmc(bayesescom3)
HPDinterval(x3,0.9)
xall3<-do.call(rbind,x3)
HPDinterval(mcmc(xall3),0.9)
x4<-as.mcmc(bayesescom4)
HPDinterval(x4,0.9)
xall4<-do.call(rbind,x4)
HPDinterval(mcmc(xall4),0.9)
x5<-as.mcmc(bayesescom5)
HPDinterval(x5,0.9)
xall5<-do.call(rbind,x5)
HPDinterval(mcmc(xall5),0.9)
x6<-as.mcmc(bayesescom6)
HPDinterval(x6,0.9)
xall6<-do.call(rbind,x6)
HPDinterval(mcmc(xall6),0.9)
x7<-as.mcmc(bayesescom7)
HPDinterval(x7,0.9)
xall7<-do.call(rbind,x7)
HPDinterval(mcmc(xall7),0.9)
x8<-as.mcmc(bayesescom8)
HPDinterval(x8,0.9)
xall8<-do.call(rbind,x8)
HPDinterval(mcmc(xall8),0.9)
x9<-as.mcmc(bayesescom9)
HPDinterval(x9,0.9)
xall9<-do.call(rbind,x9)
HPDinterval(mcmc(xall9),0.9)
x10<-as.mcmc(bayesescom10)
HPDinterval(x10,0.9)
xall10<-do.call(rbind,x10)
HPDinterval(mcmc(xall10),0.9)
x11<-as.mcmc(bayesescom11)
HPDinterval(x11,0.9)
xall11<-do.call(rbind,x11)
HPDinterval(mcmc(xall11),0.9)
x12<-as.mcmc(bayesescom12)
HPDinterval(x12,0.9)
xall12<-do.call(rbind,x12)
HPDinterval(mcmc(xall12),0.9)
x13<-as.mcmc(bayesescom13)
HPDinterval(x13,0.9)
xall13<-do.call(rbind,x13)
HPDinterval(mcmc(xall13),0.9)
x14<-as.mcmc(bayesescom14)
HPDinterval(x14,0.9)
xall14<-do.call(rbind,x14)
HPDinterval(mcmc(xall14),0.9)
x15<-as.mcmc(bayesescom15)
HPDinterval(x15,0.9)
xall15<-do.call(rbind,x15)
HPDinterval(mcmc(xall15),0.9)
x16<-as.mcmc(bayesescom16)
HPDinterval(x16,0.9)
xall16<-do.call(rbind,x16)
HPDinterval(mcmc(xall16),0.9)
x17<-as.mcmc(bayesescom17)
HPDinterval(x17,0.9)
xall17<-do.call(rbind,x17)
HPDinterval(mcmc(xall17),0.9)
x18<-as.mcmc(bayesescom18)
HPDinterval(x18,0.9)
xall18<-do.call(rbind,x18)
HPDinterval(mcmc(xall18),0.9)
x19<-as.mcmc(bayesescom19)
HPDinterval(x19,0.9)
xall19<-do.call(rbind,x19)
HPDinterval(mcmc(xall19),0.9)
x20<-as.mcmc(bayesescom20)
HPDinterval(x20,0.9)
xall20<-do.call(rbind,x20)
HPDinterval(mcmc(xall20),0.9)
x21<-as.mcmc(bayesescom21)
HPDinterval(x21,0.9)
xall21<-do.call(rbind,x21)
HPDinterval(mcmc(xall21),0.9)
x22<-as.mcmc(bayesescom22)
HPDinterval(x22,0.9)
xall22<-do.call(rbind,x22)
HPDinterval(mcmc(xall22),0.9)
x23<-as.mcmc(bayesescom23)
HPDinterval(x23,0.9)
xall23<-do.call(rbind,x23)
HPDinterval(mcmc(xall23),0.9)
x24<-as.mcmc(bayesescom24)
HPDinterval(x24,0.9)
xall24<-do.call(rbind,x24)
HPDinterval(mcmc(xall24),0.9)
x25<-as.mcmc(bayesescom25)
HPDinterval(x25,0.9)
xall25<-do.call(rbind,x25)
HPDinterval(mcmc(xall25),0.9)
x26<-as.mcmc(bayesescom26)
HPDinterval(x26,0.9)
xall26<-do.call(rbind,x26)
HPDinterval(mcmc(xall26),0.9)
x27<-as.mcmc(bayesescom27)
HPDinterval(x27,0.9)
xall27<-do.call(rbind,x27)
HPDinterval(mcmc(xall27),0.9)
x28<-as.mcmc(bayesescom28)
HPDinterval(x28,0.9)
xall28<-do.call(rbind,x28)
HPDinterval(mcmc(xall28),0.9)
x29<-as.mcmc(bayesescom29)
HPDinterval(x29,0.9)
xall29<-do.call(rbind,x29)
HPDinterval(mcmc(xall29),0.9)
x30<-as.mcmc(bayesescom30)
HPDinterval(x30,0.9)
xall30<-do.call(rbind,x30)
HPDinterval(mcmc(xall30),0.9)
x31<-as.mcmc(bayesescom31)
HPDinterval(x31,0.9)
xall31<-do.call(rbind,x31)
HPDinterval(mcmc(xall31),0.9)
x32<-as.mcmc(bayesescom32)
HPDinterval(x32,0.9)
xall32<-do.call(rbind,x32)
HPDinterval(mcmc(xall32),0.9)
x33<-as.mcmc(bayesescom33)
HPDinterval(x33,0.9)
xall33<-do.call(rbind,x33)
HPDinterval(mcmc(xall33),0.9)
x34<-as.mcmc(bayesescom34)
HPDinterval(x34,0.9)
xall34<-do.call(rbind,x34)
HPDinterval(mcmc(xall34),0.9)
x35<-as.mcmc(bayesescom35)
HPDinterval(x35,0.9)
xall35<-do.call(rbind,x35)
HPDinterval(mcmc(xall35),0.9)
x36<-as.mcmc(bayesescom36)
HPDinterval(x36,0.9)
xall36<-do.call(rbind,x36)
HPDinterval(mcmc(xall36),0.9)
x37<-as.mcmc(bayesescom37)
HPDinterval(x37,0.9)
xall37<-do.call(rbind,x37)
HPDinterval(mcmc(xall37),0.9)
x38<-as.mcmc(bayesescom38)
HPDinterval(x38,0.9)
xall38<-do.call(rbind,x38)
HPDinterval(mcmc(xall38),0.9)
x39<-as.mcmc(bayesescom39)
HPDinterval(x39,0.9)
xall39<-do.call(rbind,x39)
HPDinterval(mcmc(xall39),0.9)
x40<-as.mcmc(bayesescom40)
HPDinterval(x40,0.9)
xall40<-do.call(rbind,x40)
HPDinterval(mcmc(xall40),0.9)
x41<-as.mcmc(bayesescom41)
HPDinterval(x41,0.9)
xall41<-do.call(rbind,x41)
HPDinterval(mcmc(xall41),0.9)
x42<-as.mcmc(bayesescom42)
HPDinterval(x42,0.9)
xall42<-do.call(rbind,x42)
HPDinterval(mcmc(xall42),0.9)
x43<-as.mcmc(bayesescom43)
HPDinterval(x43,0.9)
xall43<-do.call(rbind,x43)
HPDinterval(mcmc(xall43),0.9)
x44<-as.mcmc(bayesescom44)
HPDinterval(x44,0.9)
xall44<-do.call(rbind,x44)
HPDinterval(mcmc(xall44),0.9)
x45<-as.mcmc(bayesescom45)
HPDinterval(x45,0.9)
xall45<-do.call(rbind,x45)
HPDinterval(mcmc(xall45),0.9)
x46<-as.mcmc(bayesescom46)
HPDinterval(x46,0.9)
xall46<-do.call(rbind,x46)
HPDinterval(mcmc(xall46),0.9)
x47<-as.mcmc(bayesescom47)
HPDinterval(x47,0.9)
xall47<-do.call(rbind,x47)
HPDinterval(mcmc(xall47),0.9)
x48<-as.mcmc(bayesescom48)
HPDinterval(x48,0.9)
xall48<-do.call(rbind,x48)
HPDinterval(mcmc(xall48),0.9)
x49<-as.mcmc(bayesescom49)
HPDinterval(x49,0.9)
xall49<-do.call(rbind,x49)
HPDinterval(mcmc(xall49),0.9)
x50<-as.mcmc(bayesescom50)
HPDinterval(x50,0.9)
xall50<-do.call(rbind,x50)
HPDinterval(mcmc(xall50),0.9)
x51<-as.mcmc(bayesescom51)
HPDinterval(x51,0.9)
xall51<-do.call(rbind,x51)
HPDinterval(mcmc(xall51),0.9)
x52<-as.mcmc(bayesescom52)
HPDinterval(x52,0.9)
xall52<-do.call(rbind,x52)
HPDinterval(mcmc(xall52),0.9)

# Line 8 for Expedited Coethnicity #

plot(NA, xlim=c(-5,5), ylim=c(1,10), xlab="Expedited_Coethnicity 90% Highest Posterior Density Intervals", ylab="Prior Variance")
segments(0,0,0,10,col=c("red"))
segments(HPDinterval(mcmc(xall3),.9)[8,1],10,HPDinterval(mcmc(xall3),.9)[8,2],10,col=c("grey0"))
segments(HPDinterval(mcmc(xall4),.9)[8,1],9.8,HPDinterval(mcmc(xall4),.9)[8,2],9.8,col=c("grey2"))
segments(HPDinterval(mcmc(xall5),.9)[8,1],9.6,HPDinterval(mcmc(xall5),.9)[8,2],9.6,col=c("grey 4"))
segments(HPDinterval(mcmc(xall6),.9)[8,1],9.4,HPDinterval(mcmc(xall6),.9)[8,2],9.4,col=c("grey6"))
segments(HPDinterval(mcmc(xall7),.9)[8,1],9.2,HPDinterval(mcmc(xall7),.9)[8,2],9.2,col=c("grey8"))
segments(HPDinterval(mcmc(xall8),.9)[8,1],9,HPDinterval(mcmc(xall8),.9)[8,2],9,col=c("grey10"))
segments(HPDinterval(mcmc(xall9),.9)[8,1],8.8,HPDinterval(mcmc(xall9),.9)[8,2],8.8,col=c("grey12"))
segments(HPDinterval(mcmc(xall10),.9)[8,1],8.6,HPDinterval(mcmc(xall10),.9)[8,2],8.6,col=c("grey14"))
segments(HPDinterval(mcmc(xall11),.9)[8,1],8.40000000000001,HPDinterval(mcmc(xall11),.9)[8,2],8.40000000000001,col=c("grey16"))
segments(HPDinterval(mcmc(xall12),.9)[8,1],8.20000000000001,HPDinterval(mcmc(xall12),.9)[8,2],8.20000000000001,col=c("grey18"))
segments(HPDinterval(mcmc(xall13),.9)[8,1],8.00000000000001,HPDinterval(mcmc(xall13),.9)[8,2],8.00000000000001,col=c("grey20"))
segments(HPDinterval(mcmc(xall14),.9)[8,1],7.80000000000001,HPDinterval(mcmc(xall14),.9)[8,2],7.80000000000001,col=c("grey22"))
segments(HPDinterval(mcmc(xall15),.9)[8,1],7.60000000000001,HPDinterval(mcmc(xall15),.9)[8,2],7.60000000000001,col=c("grey24"))
segments(HPDinterval(mcmc(xall16),.9)[8,1],7.40000000000001,HPDinterval(mcmc(xall16),.9)[8,2],7.40000000000001,col=c("grey26"))
segments(HPDinterval(mcmc(xall17),.9)[8,1],7.20000000000001,HPDinterval(mcmc(xall17),.9)[8,2],7.20000000000001,col=c("grey28"))
segments(HPDinterval(mcmc(xall18),.9)[8,1],7.00000000000001,HPDinterval(mcmc(xall18),.9)[8,2],7.00000000000001,col=c("grey30"))
segments(HPDinterval(mcmc(xall19),.9)[8,1],6.80000000000001,HPDinterval(mcmc(xall19),.9)[8,2],6.80000000000001,col=c("grey32"))
segments(HPDinterval(mcmc(xall20),.9)[8,1],6.60000000000001,HPDinterval(mcmc(xall20),.9)[8,2],6.60000000000001,col=c("grey34"))
segments(HPDinterval(mcmc(xall21),.9)[8,1],6.40000000000001,HPDinterval(mcmc(xall21),.9)[8,2],6.40000000000001,col=c("grey36"))
segments(HPDinterval(mcmc(xall22),.9)[8,1],6.20000000000001,HPDinterval(mcmc(xall22),.9)[8,2],6.20000000000001,col=c("grey38"))
segments(HPDinterval(mcmc(xall23),.9)[8,1],6.00000000000001,HPDinterval(mcmc(xall23),.9)[8,2],6.00000000000001,col=c("grey40"))
segments(HPDinterval(mcmc(xall24),.9)[8,1],5.80000000000001,HPDinterval(mcmc(xall24),.9)[8,2],5.80000000000001,col=c("grey42"))
segments(HPDinterval(mcmc(xall25),.9)[8,1],5.60000000000002,HPDinterval(mcmc(xall25),.9)[8,2],5.60000000000002,col=c("grey44"))
segments(HPDinterval(mcmc(xall26),.9)[8,1],5.40000000000002,HPDinterval(mcmc(xall26),.9)[8,2],5.40000000000002,col=c("grey46"))
segments(HPDinterval(mcmc(xall27),.9)[8,1],5.20000000000002,HPDinterval(mcmc(xall27),.9)[8,2],5.20000000000002,col=c("grey48"))
segments(HPDinterval(mcmc(xall28),.9)[8,1],5.00000000000002,HPDinterval(mcmc(xall28),.9)[8,2],5.00000000000002,col=c("grey50"))
segments(HPDinterval(mcmc(xall29),.9)[8,1],4.80000000000002,HPDinterval(mcmc(xall29),.9)[8,2],4.80000000000002,col=c("grey52"))
segments(HPDinterval(mcmc(xall30),.9)[8,1],4.60000000000002,HPDinterval(mcmc(xall30),.9)[8,2],4.60000000000002,col=c("grey54"))
segments(HPDinterval(mcmc(xall31),.9)[8,1],4.40000000000002,HPDinterval(mcmc(xall31),.9)[8,2],4.40000000000002,col=c("grey56"))
segments(HPDinterval(mcmc(xall32),.9)[8,1],4.20000000000002,HPDinterval(mcmc(xall32),.9)[8,2],4.20000000000002,col=c("grey58"))
segments(HPDinterval(mcmc(xall33),.9)[8,1],4.00000000000002,HPDinterval(mcmc(xall33),.9)[8,2],4.00000000000002,col=c("grey60"))
segments(HPDinterval(mcmc(xall34),.9)[8,1],3.80000000000002,HPDinterval(mcmc(xall34),.9)[8,2],3.80000000000002,col=c("grey62"))
segments(HPDinterval(mcmc(xall35),.9)[8,1],3.60000000000002,HPDinterval(mcmc(xall35),.9)[8,2],3.60000000000002,col=c("grey64"))
segments(HPDinterval(mcmc(xall36),.9)[8,1],3.40000000000002,HPDinterval(mcmc(xall36),.9)[8,2],3.40000000000002,col=c("grey66"))
segments(HPDinterval(mcmc(xall37),.9)[8,1],3.20000000000002,HPDinterval(mcmc(xall37),.9)[8,2],3.20000000000002,col=c("grey68"))
segments(HPDinterval(mcmc(xall38),.9)[8,1],3.00000000000002,HPDinterval(mcmc(xall38),.9)[8,2],3.00000000000002,col=c("grey70"))
segments(HPDinterval(mcmc(xall39),.9)[8,1],2.80000000000003,HPDinterval(mcmc(xall39),.9)[8,2],2.80000000000003,col=c("grey72"))
segments(HPDinterval(mcmc(xall40),.9)[8,1],2.60000000000003,HPDinterval(mcmc(xall40),.9)[8,2],2.60000000000003,col=c("grey74"))
segments(HPDinterval(mcmc(xall41),.9)[8,1],2.40000000000003,HPDinterval(mcmc(xall41),.9)[8,2],2.40000000000003,col=c("grey76"))
segments(HPDinterval(mcmc(xall42),.9)[8,1],2.20000000000003,HPDinterval(mcmc(xall42),.9)[8,2],2.20000000000003,col=c("grey78"))
segments(HPDinterval(mcmc(xall43),.9)[8,1],2.00000000000003,HPDinterval(mcmc(xall43),.9)[8,2],2.00000000000003,col=c("grey80"))
segments(HPDinterval(mcmc(xall44),.9)[8,1],1.80000000000003,HPDinterval(mcmc(xall44),.9)[8,2],1.80000000000003,col=c("grey82"))
segments(HPDinterval(mcmc(xall45),.9)[8,1],1.60000000000003,HPDinterval(mcmc(xall45),.9)[8,2],1.60000000000003,col=c("grey84"))
segments(HPDinterval(mcmc(xall46),.9)[8,1],1.40000000000003,HPDinterval(mcmc(xall46),.9)[8,2],1.40000000000003,col=c("grey86"))
segments(HPDinterval(mcmc(xall47),.9)[8,1],1.20000000000003,HPDinterval(mcmc(xall47),.9)[8,2],1.20000000000003,col=c("grey88"))
segments(HPDinterval(mcmc(xall48),.9)[8,1],1.00000000000003,HPDinterval(mcmc(xall48),.9)[8,2],1.00000000000003,col=c("grey90"))
segments(HPDinterval(mcmc(xall49),.9)[8,1],0.800000000000029,HPDinterval(mcmc(xall49),.9)[8,2],0.800000000000029,col=c("grey92"))
segments(HPDinterval(mcmc(xall50),.9)[8,1],0.60000000000003,HPDinterval(mcmc(xall50),.9)[8,2],0.60000000000003,col=c("grey94"))
segments(HPDinterval(mcmc(xall51),.9)[8,1],0.400000000000031,HPDinterval(mcmc(xall51),.9)[8,2],0.400000000000031,col=c("grey96"))
segments(HPDinterval(mcmc(xall52),.9)[8,1],0.200000000000029,HPDinterval(mcmc(xall52),.9)[8,2],0.200000000000029,col=c("grey98"))

# Line 3 for Bribery Coethnicity #

plot(NA, xlim=c(-5,5), ylim=c(1,10), xlab="Bribery_Coethnicity 90% Highest Posterior Density Intervals", ylab="Prior Variance")
segments(0,0,0,10,col=c("red"))
segments(HPDinterval(mcmc(xall3),.9)[3,1],10,HPDinterval(mcmc(xall3),.9)[3,2],10,col=c("grey0"))
segments(HPDinterval(mcmc(xall4),.9)[3,1],9.8,HPDinterval(mcmc(xall4),.9)[3,2],9.8,col=c("grey2"))
segments(HPDinterval(mcmc(xall5),.9)[3,1],9.6,HPDinterval(mcmc(xall5),.9)[3,2],9.6,col=c("grey 4"))
segments(HPDinterval(mcmc(xall6),.9)[3,1],9.4,HPDinterval(mcmc(xall6),.9)[3,2],9.4,col=c("grey6"))
segments(HPDinterval(mcmc(xall7),.9)[3,1],9.2,HPDinterval(mcmc(xall7),.9)[3,2],9.2,col=c("grey8"))
segments(HPDinterval(mcmc(xall8),.9)[3,1],9,HPDinterval(mcmc(xall8),.9)[3,2],9,col=c("grey10"))
segments(HPDinterval(mcmc(xall9),.9)[3,1],8.8,HPDinterval(mcmc(xall9),.9)[3,2],8.8,col=c("grey12"))
segments(HPDinterval(mcmc(xall10),.9)[3,1],8.6,HPDinterval(mcmc(xall10),.9)[3,2],8.6,col=c("grey14"))
segments(HPDinterval(mcmc(xall11),.9)[3,1],8.40000000000001,HPDinterval(mcmc(xall11),.9)[3,2],8.40000000000001,col=c("grey16"))
segments(HPDinterval(mcmc(xall12),.9)[3,1],8.20000000000001,HPDinterval(mcmc(xall12),.9)[3,2],8.20000000000001,col=c("grey18"))
segments(HPDinterval(mcmc(xall13),.9)[3,1],8.00000000000001,HPDinterval(mcmc(xall13),.9)[3,2],8.00000000000001,col=c("grey20"))
segments(HPDinterval(mcmc(xall14),.9)[3,1],7.80000000000001,HPDinterval(mcmc(xall14),.9)[3,2],7.80000000000001,col=c("grey22"))
segments(HPDinterval(mcmc(xall15),.9)[3,1],7.60000000000001,HPDinterval(mcmc(xall15),.9)[3,2],7.60000000000001,col=c("grey24"))
segments(HPDinterval(mcmc(xall16),.9)[3,1],7.40000000000001,HPDinterval(mcmc(xall16),.9)[3,2],7.40000000000001,col=c("grey26"))
segments(HPDinterval(mcmc(xall17),.9)[3,1],7.20000000000001,HPDinterval(mcmc(xall17),.9)[3,2],7.20000000000001,col=c("grey28"))
segments(HPDinterval(mcmc(xall18),.9)[3,1],7.00000000000001,HPDinterval(mcmc(xall18),.9)[3,2],7.00000000000001,col=c("grey30"))
segments(HPDinterval(mcmc(xall19),.9)[3,1],6.80000000000001,HPDinterval(mcmc(xall19),.9)[3,2],6.80000000000001,col=c("grey32"))
segments(HPDinterval(mcmc(xall20),.9)[3,1],6.60000000000001,HPDinterval(mcmc(xall20),.9)[3,2],6.60000000000001,col=c("grey34"))
segments(HPDinterval(mcmc(xall21),.9)[3,1],6.40000000000001,HPDinterval(mcmc(xall21),.9)[3,2],6.40000000000001,col=c("grey36"))
segments(HPDinterval(mcmc(xall22),.9)[3,1],6.20000000000001,HPDinterval(mcmc(xall22),.9)[3,2],6.20000000000001,col=c("grey38"))
segments(HPDinterval(mcmc(xall23),.9)[3,1],6.00000000000001,HPDinterval(mcmc(xall23),.9)[3,2],6.00000000000001,col=c("grey40"))
segments(HPDinterval(mcmc(xall24),.9)[3,1],5.80000000000001,HPDinterval(mcmc(xall24),.9)[3,2],5.80000000000001,col=c("grey42"))
segments(HPDinterval(mcmc(xall25),.9)[3,1],5.60000000000002,HPDinterval(mcmc(xall25),.9)[3,2],5.60000000000002,col=c("grey44"))
segments(HPDinterval(mcmc(xall26),.9)[3,1],5.40000000000002,HPDinterval(mcmc(xall26),.9)[3,2],5.40000000000002,col=c("grey46"))
segments(HPDinterval(mcmc(xall27),.9)[3,1],5.20000000000002,HPDinterval(mcmc(xall27),.9)[3,2],5.20000000000002,col=c("grey48"))
segments(HPDinterval(mcmc(xall28),.9)[3,1],5.00000000000002,HPDinterval(mcmc(xall28),.9)[3,2],5.00000000000002,col=c("grey50"))
segments(HPDinterval(mcmc(xall29),.9)[3,1],4.80000000000002,HPDinterval(mcmc(xall29),.9)[3,2],4.80000000000002,col=c("grey52"))
segments(HPDinterval(mcmc(xall30),.9)[3,1],4.60000000000002,HPDinterval(mcmc(xall30),.9)[3,2],4.60000000000002,col=c("grey54"))
segments(HPDinterval(mcmc(xall31),.9)[3,1],4.40000000000002,HPDinterval(mcmc(xall31),.9)[3,2],4.40000000000002,col=c("grey56"))
segments(HPDinterval(mcmc(xall32),.9)[3,1],4.20000000000002,HPDinterval(mcmc(xall32),.9)[3,2],4.20000000000002,col=c("grey58"))
segments(HPDinterval(mcmc(xall33),.9)[3,1],4.00000000000002,HPDinterval(mcmc(xall33),.9)[3,2],4.00000000000002,col=c("grey60"))
segments(HPDinterval(mcmc(xall34),.9)[3,1],3.80000000000002,HPDinterval(mcmc(xall34),.9)[3,2],3.80000000000002,col=c("grey62"))
segments(HPDinterval(mcmc(xall35),.9)[3,1],3.60000000000002,HPDinterval(mcmc(xall35),.9)[3,2],3.60000000000002,col=c("grey64"))
segments(HPDinterval(mcmc(xall36),.9)[3,1],3.40000000000002,HPDinterval(mcmc(xall36),.9)[3,2],3.40000000000002,col=c("grey66"))
segments(HPDinterval(mcmc(xall37),.9)[3,1],3.20000000000002,HPDinterval(mcmc(xall37),.9)[3,2],3.20000000000002,col=c("grey68"))
segments(HPDinterval(mcmc(xall38),.9)[3,1],3.00000000000002,HPDinterval(mcmc(xall38),.9)[3,2],3.00000000000002,col=c("grey70"))
segments(HPDinterval(mcmc(xall39),.9)[3,1],2.80000000000003,HPDinterval(mcmc(xall39),.9)[3,2],2.80000000000003,col=c("grey72"))
segments(HPDinterval(mcmc(xall40),.9)[3,1],2.60000000000003,HPDinterval(mcmc(xall40),.9)[3,2],2.60000000000003,col=c("grey74"))
segments(HPDinterval(mcmc(xall41),.9)[3,1],2.40000000000003,HPDinterval(mcmc(xall41),.9)[3,2],2.40000000000003,col=c("grey76"))
segments(HPDinterval(mcmc(xall42),.9)[3,1],2.20000000000003,HPDinterval(mcmc(xall42),.9)[3,2],2.20000000000003,col=c("grey78"))
segments(HPDinterval(mcmc(xall43),.9)[3,1],2.00000000000003,HPDinterval(mcmc(xall43),.9)[3,2],2.00000000000003,col=c("grey80"))
segments(HPDinterval(mcmc(xall44),.9)[3,1],1.80000000000003,HPDinterval(mcmc(xall44),.9)[3,2],1.80000000000003,col=c("grey82"))
segments(HPDinterval(mcmc(xall45),.9)[3,1],1.60000000000003,HPDinterval(mcmc(xall45),.9)[3,2],1.60000000000003,col=c("grey84"))
segments(HPDinterval(mcmc(xall46),.9)[3,1],1.40000000000003,HPDinterval(mcmc(xall46),.9)[3,2],1.40000000000003,col=c("grey86"))
segments(HPDinterval(mcmc(xall47),.9)[3,1],1.20000000000003,HPDinterval(mcmc(xall47),.9)[3,2],1.20000000000003,col=c("grey88"))
segments(HPDinterval(mcmc(xall48),.9)[3,1],1.00000000000003,HPDinterval(mcmc(xall48),.9)[3,2],1.00000000000003,col=c("grey90"))
segments(HPDinterval(mcmc(xall49),.9)[3,1],0.800000000000029,HPDinterval(mcmc(xall49),.9)[3,2],0.800000000000029,col=c("grey92"))
segments(HPDinterval(mcmc(xall50),.9)[3,1],0.60000000000003,HPDinterval(mcmc(xall50),.9)[3,2],0.60000000000003,col=c("grey94"))
segments(HPDinterval(mcmc(xall51),.9)[3,1],0.400000000000031,HPDinterval(mcmc(xall51),.9)[3,2],0.400000000000031,col=c("grey96"))
segments(HPDinterval(mcmc(xall52),.9)[3,1],0.200000000000029,HPDinterval(mcmc(xall52),.9)[3,2],0.200000000000029,col=c("grey98"))